[国家集训队]Crash的数字表格 / JZPTAB

可以知道,题目中所求的是

i=1nj=1mlcm(i,j)=i=1nj=1mijgcd(i,j)\sum\limits_{i=1}^n\sum\limits_{j=1}^mlcm(i,j)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{ij}{gcd(i,j)}

枚举所有约数

i=1nj=1mijgcd(i,j)=d=1min(n,m)i=1nj=1mijd[gcd(i,j)=d]=d=1min(n,m)di=1ndj=1mdij[gcd(i,j)=1]\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{ij}{gcd(i,j)}=\sum\limits_{d=1}^{min(n,m)}\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{ij}{d}[gcd(i,j)=d]\\=\sum\limits_{d=1}^{min(n,m)}d\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij[gcd(i,j)=1]

sum(n,m)=i=1nj=1mij[gcd(i,j)=1]f(k)=i=1nj=1mij[gcd(i,j)=k],F(k)=i=1nj=1mij[kgcd(i,j)]sum(n,m)=\sum\limits_{i=1}^n\sum\limits_{j=1}^mij[gcd(i,j)=1]\\f(k)=\sum\limits_{i=1}^n\sum\limits_{j=1}^mij[gcd(i,j)=k],F(k)=\sum\limits_{i=1}^n\sum\limits_{j=1}^mij[k|gcd(i,j)]

F(k)=kdf(d)f(k)=kdμ(dk)F(d)=kdμ(dk)i=1nj=1mij[dgcd(i,j)]=kdμ(dk)d2i=1ndj=1mdij\Rightarrow F(k)=\sum\limits_{k|d}f(d)\Rightarrow f(k)=\sum\limits_{k|d}\mu(\lfloor\frac{d}{k}\rfloor)F(d)\\=\sum\limits_{k|d}\mu(\lfloor\frac{d}{k}\rfloor)\sum\limits_{i=1}^n\sum\limits_{j=1}^mij[d|gcd(i,j)]=\sum\limits_{k|d}\mu(\lfloor\frac{d}{k}\rfloor)d^2\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij

g(i,j)=i×(i+1)2×j×(j+1)2g(i,j)=\frac{i\times(i+1)}{2}\times\frac{j\times(j+1)}{2}

f(k)=kdμ(dk)d2g(nd,md)\Rightarrow f(k)=\sum\limits_{k|d}\mu(\lfloor\frac{d}{k}\rfloor)d^2g(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)

sum(n,m)=f(1)=d=1min(n,m)μ(d)d2g(nd,md)\Rightarrow sum(n,m)=f(1)=\sum\limits_{d=1}^{min(n,m)}\mu(d)d^2g(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)

ans=d=1min(n,m)dsum(nd,md)ans=\sum\limits_{d=1}^{min(n,m)}dsum(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)

推到这就很明显可以用整除分块了,时间复杂度 O(n)O(n),注意取模和long long

代码

#include<iostream>
#include<cstdio>
#define MOD (20101009)
#define g(x, y) ((((ll)(x + 1) * (x) / 2 % MOD) * ((ll)(y + 1) * (y) / 2 % MOD)) % MOD)
using namespace std;
const int N = 1e7 + 5;
typedef long long ll;

bool st[N];
int prime[N], mu[N];
ll s[N];

void init()
{
int tot = 0;
mu[1] = 1;
for (int i = 2; i < N; i++)
{
if (!st[i])
prime[++tot] = i, mu[i] = -1;
for (int j = 1; i * prime[j] < N; j++)
{
st[i * prime[j]] = true;
if (i % prime[j] == 0)
break;
mu[i * prime[j]] = -mu[i];
}
}
for (int i = 1; i < N; i++)
s[i] = ((ll)s[i - 1] + (ll)mu[i] * i * i + MOD) % MOD;
}

ll sum(int n, int m)
{
int k = min(n, m);
ll res = 0;
for (int l = 1, r; l <= k; l = r + 1)
{
r = min(k, min(n / (n / l), m / (m / l)));
res = (res + (ll)(s[r] - s[l - 1] + MOD) * g(n / l, m / l) % MOD) % MOD;
}
return res;
}

int main()
{
init();
int n, m;
scanf("%d%d", &n, &m);
int k = min(n, m);
ll res = 0;
for (int l = 1, r; l <= k; l = r + 1)
{
r = min(k, min(n / (n / l), m / (m / l)));
res += (ll)(l + r) * (r - l + 1) / 2 * sum(n / l, m / l);
res %= MOD;
}
printf("%lld", res);
return 0;
}